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ABSTRACT 

By stacking an ensemble of strong lensing clusters, we demonstrate the feasibility of 
placing constraints on the dark energy equation of state. This is achieved by using mul- 
tiple images of sources at two or more distinct redshift planes. The sample of smooth 
clusters in our simulations is based on observations of massive clusters and the distri- 
bution of background galaxies is constructed using the Hubble Deep Field. Our source 
distribution reproduces the observed redshift distribution of multiply imaged sources 
in Abell 1689. The cosmology recovery depends on the number of image families with 
known spectroscopic redshifts and the number of stacked clusters. Our simulations 
suggest that constraints comparable to those derived from other competing estab- 
lished techniques on a constant dark energy equation of state can be obtained using 
10 to 40 clusters with 5 or more families of multiple images. We have also studied the 
observational errors in the image redshifts and positions. We find that spectroscopic 
redshifts and high resolution Hubble Space Telescope images are required to eliminate 
confidence contour relaxation relative to the ideal case in our simulations. This sug- 
gests that the dark energy equation of state, and other cosmological parameters, can 
be constrained with existing Hubble Space Telescope images of lensing clusters coupled 
with dedicated ground-based arc spectroscopy. 
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1 INTRODUCTION 

Current efforts in observational cosmology are directed 
towards characterizing the energy content of the Uni- 
verse. The discovery of the accelerating expansion of 
the Universe inferred fr o m th e Hubble diagram for 
SN type la jRiess et al l Il998l; [Perlmutter et al.1 1 19991 ; 



iTonrv et all 120031 ; iRiess et all |2004 120071 ) combined with 
constraints from Cosmic Microwave Background Radia- 
tion from th e Wilk i nson Microwave Anisotropy Probe 
(|Spergel et all l200l 120071 ; iHinshaw et al l 120091) from 



cosmic shear observations llBacon. Refregier fc Ellis! 200C 



Kaiseij|2000l;lvan Waerbeke et al.ll2000l : IWittman et allboOC 
Semboloni et al.ll2006l ), cluster baryon fractions \ Allen et al 
20041 ') . galaxy surveys (|Efstathiou et al" 1 120021 ; ISeliak et al 
2005) suggests that ~ 70 per cent of the energy in the Uni- 



verse is in the form of dark energy. Dark energy is best 
described by the equation of state that relates its pressure, 
P, to its energy density, p, via P — wxp, with wx = — 1 
corresponding to the case of the cosmological constant. 

Massive foreground clusters often produce many sets 
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of highly distorted arcs and multiple images of background 
galaxies ( Blan dford fc Naravan|[l992l ). The lensing effect is 
determined by the mass distribution, the precise configura- 
tion of the lens and source with respect to the observer, and 
the angular diameter distances between the lens, source and 
observer. When at least two sources at distinct redshift are 
strongly lensed in a cluster, the relative positions of the arcs 
depend on the ratios of angular diameter distances. In prin- 
ciple therefore, cosmological parameters can be constrained 
from the lensing configuration. 

Constraining the geometry and matter content of the 
Universe using multiple sets of arcs in cluster lenses has 
been e x plored by Paczynski & Gorski (1981); Link fc Pierce! 
(1 1998!) ; ICooravi (1 19991); iGolse. Kneib fc Soucaill (|2002T ): 
Serend I2OO2I); ISoucail. Kneib fc Golsel i|2004l ). Recently. 
Biesiadal (|2006l ) has suggested using multiple images 
produced by galaxy lenses instead of clusters and at- 
tempted an applic ation to the c a se of the lensing system 
HST 14176+5226. iLink fc Pierce] (Il998l ) showed that the 
cosmological sensitivity of the angular size-redshift relation 
could be exploited using sources at distinct redshifts and 
developed a methodology to simulta neously invert the lens 
and derive cosmological constraints. ISeren 3 (|2002l ) applied 
the technique to the cluster CI 0024+1654 and found sup- 
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port fo r an accelerating Universe. iGolse. Kneib. fc Soucaill 
l|2002l ) showed that the recovery of cosmological parameters 
was feasible wi th 3 sets of multiple images fo r a single sim- 
ulated cluster. ISoucail. Kneib fc Golsd (|2004h then applied 
the technique to the lensing cluster Abell 2218 and using 4 
multiple image systems at distinct redshifts, find Qm < 0.37 
and w < —0.80 assuming a flat Universe. Recently, a dou- 
ble Einstei n ring has been dis covered and applied to cos- 
mography l|Gavazzi et alJl200a ) . As the system has only two 
sources, only weak cosmography constraints were obtained. 

There have been several studies that have explored 
the cross-sections of multiple image formation using sim- 
ulated clusters f rom c osmolo gical jV-body simulations. 
iMeneehetti et all (|2005al lbh and iMacciol (|2005l ) studied the 
efficiency of numerical clusters to produce multiple images 
in different dark energy models. They found the expected 
abundance of arcs with a given length-to-width ratio and 
the lensing cross-sections depend on the equation of state 
of dark energy. These statistical studies of multiple image 
formation concluded that strong lensing offers a plausible 
method to discriminate amongst various dark energy mod- 
els. 

In this paper we study two basic aspects of the non- 
statistical approach to cosmography with cluster strong lens- 
ing (CSL), which have yet to be addressed. These are: 1) 
How many clusters must be individually marginalized, with 
the resulting \ 2 from each cluster added together to form a 
final x 2 f° r cosmology recovery (i.e. stacked), in the ideal 
case to obtain cosmological constraints comparable with 
other methods and 2) The required accuracy of the ob- 
servational parameters, i.e. the source redshift and image 
positions. By studying these two aspects separately from 
the subtleties involved with cluster mass distribution recon- 
struction, we can make definitive statements on the type of 
cluster surveys required for a first application of CSL. 

We do not model the lensing effects induced by 
the distribution of substructure along the line of sight. 
This addition al complication has been exam i ned in ear- 
lier work by Warnbsganss, Bode & Ostriker (2005) and 
IDalal. Hennawi fc Bodd l|2005h . who both used ray trac- 
ing methods through A r -body simulations. Simulating a 
representative shell of the Un i verse centered at z = 1, 
IWambsganss. Bode fc Ostrikerl (|2005l ) found that approxi- 
mately 10%-15% of the lensing mass was to be found in 
substructure along the line of sight, for about one third of 
lensing configurations with a source redshift above z 3 — 2.5. 
In particular, it was shown that substructure along the line 
of sight is limited to a finite number of lens planes, typically 
one or two additional planes, when line of sight substruc- 
ture played a significant role in lensing. They have shown 
explicitly that < 20% of the contribution to the convergence 
arises from secondary lens planes. For cluster strong lensing 
the critical factor is not the error induced in the conver- 
gence but rather the perturbations caused in the positions 
of the multiple images due to the presence of matter in ad- 
ditional lens planes. A detailed and exhaustive set of cosmo- 
logical simulations that will enable the matter distribution 
of a typical line of sight to be averaged over many repre- 
sentative lines of sight is needed t o accu rately quantify this 
effect. In lDalal. Hennawi fc Bode] (120051 ) it was shown that 
systematic errors in the cosmology recovery could result, if 
the line of sight substructure was not taken into account. 



The outline of this paper is as follows: we begin with an 
introduction to the basics of strong lensing and dark energy 
in Section 2. Then in Section 3 we discuss image families, and 
define the cosmological family ratio, which CSL constrains. 
The methodology used to generate smooth simulated clus- 
ters is described in Section 4 and Section 5, along with the 
\ 2 minimization scheme used in the cosmology recovery. We 
present the results of recovering the input cosmology in Sec- 
tion 6, and study the dependence of the recovery on lens 
parameters in Section 7. We then study the allowed errors 
in the source redshift and image positions in Section 8. This 
is followed by a discussion of the simplifying assumptions 
made in this work in Section 9. We conclude with a sum- 
mary of our technique, and the implications for future survey 
strategies. 



2 STRONG GRAVITATIONAL LENSING AND 
DARK ENERGY 

The gravitational lens equation is a mapping from the source 
plane at redshift z s to the image plane (or lens plane) at 
redshift z\ and is given by, 



= (3 + a(0,£;M) 



(1) 



where 9 is the angular position of the image, (5 is the un- 
observed source position, a. is the deflection angle which 
depends on £ the reduced angular diameter distance, M the 
mass distribution of the lens and 6. The cosmological de- 
pendence is contained in £, which we call the reduced an- 
gular diameter defined by: £ = D(0, z\)D(z\, z B )/D(Q, z s ) = 
D \ D\ s / D os , where D{z\,Z2) is the angular diameter dis- 
tance between redshift z\ and zi, and D a \, D\ s and D oa 
are the observer-lens, lens-source and observer-source an- 
gular diameter distances respectively. In a flat homogeneous 
Friedmann-Robertson- Walker cosmology the angular diam- 
eter distance is given by, 
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where c is the speed of light; Ho = 100 ft km s~ is the 
present day value of the Hubble constant, and Qm and 
f2x are the fractional contributions of the total matter and 
dark energy respectively, in units of the critical density and 
wx is the dark energy equation of state. Note that we do 
not consider evolving dark energy models in our simula- 
tions, since there is currently no observational evidence for 
such a time evolution in the dark energy equatio n of state 
llSahni fc Starobinskvll2006l : iKomatsu et alj|2009h . 

To further explore the dependence of the lens equation 
on cosmological parameters, we note that the deflection an- 
gle can be written in terms of the gradient of the reduced 
gravitational potential ipm(6,£)'- V<£Af(#,£) = cx(0,£;M). 
The reduced gravitational potential can be related to the 
projected surface mass density along the line of sight, 
Ej\/(0), through the Poisson equation, where the projected 
gravitational potential (/)m(9) is defined by 



(3) 



In the thin lens approximation the projected gravitational 
potential is a solution of: V 2 ^m(#,£) = 8nG £ Em(0)/c 2 . 
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The total deflection (G — (3) is therefore determined both 
by the cosmological parameters (via D \Di a / D os ) and the 
projected surface mass density of the intervening lens. For 
a given configuration of the source, lens and observer, the 
critical value of the surface mass density is: 



4vrG £ 



(4) 



The region where the projected surface mass density of a lens 
exceeds this critical value defines the strong lensing regime. 
Background sources that are aligned behind this region are 
likely to be multiply imaged, magnified and distorted. 

There are hence two fundamental equations which de- 
fine the gravitational lensing dependence on the angular di- 
ameter distances D a \ , D\ B , D OB . The first is the lens mapping, 
equation JTJ, which defines the configuration of the images 
produced by a lensing cluster, and the second is equation Q 
which defines the strong lensing regime. The cosmological 
dependence of strong lensing is discussed further in Sec- 
tion [3] 



3 THE COSMOLOGICAL FAMILY RATIO 
3.1 Multiple image families 

Consider a lensing configuration where at least two sources 
have been strongly lensed by a cluster. Once the multiple 
images of each source have been identified, we can write 
the source and the associated image positions with a family 
index /. The lens equation is then, 



P f = Of,i-V<pM{Of,i,t), 



(•») 



where 0f yi denotes the position of the i th image of family /. 
When there are k multiple images detected in family /, we 
can sum over all the images in equation ([5]). Substituting into 
equation (|3]) for the reduced gravitational potential allows 
the summed equation to be written in terms of the line of 
sight projected gravitational potential as follows: 



Du 



D ol D 



(6) 



Note that equation ([6]) is a non-linear function of the image 
positions, since the sum of the image positions depends on 
the projected cluster potential gradient. 

Given a set of input cosmological parameters, the lens 
redshift and a source redshift for a single family of images, 
i.e. (z\, z s , Am j f^x, wx), the angular diameter distance ratio 
is specified and constant. The cluster potential is conven- 
tionally normalized by the square of the velocity dispersion, 
which we denote by v. In equation (|6]) the reduced angu- 
lar diameter distance £ and the cluster mass normalization 
v 2 are degenerate when recovering cosmological parameters. 
Therefore, cosmological parameters cannot be disentangled 
from the cluster mass normalization with only one family 
of multiple images. Possib le solutions to remove this de - 
generacy are considered bv iGolse! Kneib. fc Soucaill (2002). 
They show an additional prior on the mass profile is re- 
quired and this must come from observations independent 
of strong lensing. For parametric models, this is typically a 
total mass determination from the detected X-ray emission 
of the cluster. 



3.2 Two families of multiple images 

Now consider the case when two strongly lensed sources are 
observed as multiple image families in the same cluster lens 
at two redshifts. Let family 1 have a total of m observed 
images and family 2 have n observed images. We require that 
at least two images from each are observed, i.e. m,n ^ 2. 
Taking the ratio of family 1 to family 2 (equation [6]) we 
obtain, 



Asl A.B2 1 Efel V0m(01 
D os l D\, 
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(7) 



Note that we now have a ratio of two factors of V</>m (#)• 
This clearly means that the degeneracy between the cosmo- 
logical parameters and the cluster potential normalization 
is broken when 2 families of multiple images are observed in 
a cluster lens. 

In equation J7), the term in braces is a function of the 
angular diameter distances only. This quantity is fundamen- 
tal to the determination of cosmological parameters through 
CSL and we define this term as the family ratio and denote 
it by H: 



H(2l, Z S 1 , Zs2 ; , fix , W X ) 



D(z h z sl ) D(0,z s2 ) 
D(0,z sl ) D(zi,z s2 ) 



(8) 



where z\ is the lens redshift, z s i and z S 2 are the two source 
redshifts, and D(zi, z 2 ) is the angular diameter distance de- 
fined in equation ([2j. The cosmological family ratio can be 
written in terms of the reduced angular diameter distance, 
S = £i/£2, where the subscripts 1 and 2 here denote the 
family index. 

We plot S in Fig. Q]for a fixed lens redshift of Z\ = 0.3 
and differing source redshifts. Since the non-evolving dark 
energy component with wx = —1 can also be interpreted 
as the cosmological constant Qa, we plot the contours of 
the cosmological family ratio in both the J^m-Oa plane in 
addition to the SIm-wx plane. Note that the contours in the 
JIm-Oa plane are not as sensitive to the lens-source config- 
uration as they are in the S7m-wx plane. 



4 SIMULATIONS OF STRONG LENSING 
CLUSTERS 

As shown in the previous section, the cosmological family 
ratio, E, depends on the source redshifts and the lens red- 
shift (see Fig. [TJ. It is therefore clear that the distribution 
of sources in magnitude and redshift as well as the clus- 
ter redshift distribution are key inputs in our simulations. 
Therefore we base our simulations of strong lensing clusters 
on high resolution observational data from HST images. We 
note that although the cluster normalization is irrelevant to 
our method for constraining the cosmological parameters, 
we still construct this quantity rigorously. This ensures the 
simulations produce lensing clusters with cross-sections sim- 
ilar to those of observed clusters. 



4.1 Source magnitude and redshift distribution 

In order to simulate the lensed images produced by a mas- 
sive cluster, it is necessary to characterize the number den- 
sity, N, of the sources behind the lens. We consider back- 



4 Gilmore & Natarajan 




Figure 1. Contours of the cosmological family ratio B(zi, z s i, Zs2', Hm> ^Xi w x) on the £7m-^A antl ^M _w x planes. The cross marks the 
input cosmology. The redshift of the lensing cluster is fixed at Z\ = 0.3. Solid contours correspond to sources at z s i = 0.7, z s 2 = 1.6 and 
dashed to sources at z a i = 2.0, z s 2 = 4.0 in both panels. 



ground point sources for which the number density is a func- 
tion of magnitude and redshift, i.e. N = N(m,z). Further- 
more, space-based positional accuracy is assumed through- 
out this paper. We must therefore draw redshifts and mag- 
nitudes from a spac e-based source distribu tion. The Hubble 



,413(8140) = 19.0 - 28.0 in intervals of 0.5 and in redshift 
from z = 0.4 — 5.0 in intervals of 0.2. These data are then 
fit with functions of the form, 



Deep Field (HD F) IWilliams et all Il99fj ) and HDF-South 
(HDF-S) WFPC2 lCasertano et all l|2000l ) are the ideal data 
sets to derive the source distribution. It is also possible to 
use the recent COSMOS data to derive the redsh i ft an d 
magnitude source distributions, see iGabasch et al. | (|2008l ). 
In particular we have checked that our source redshift dis- 
tribution, which is d erived below, is equivalent to that of 
IGabasch et ail (|200ST ). 

We employ the WFPC2 HDF and HDF-S photometric 
reds hift catalogues derived by the Stony Brook collabora- 



tion jFernandez-Soto. Lanzetta. &: Yahil 1999; Y ahata et al.l 
2000) and characterize the HDF sources by a distribution 
function. This method is chosen since we wish to analyse 
an ensemble of clusters, where each cluster lenses a unique 
source configuration. Characterizing the source distribution 
by a function is preferable in our case, since we do not want 
to bias the results by lensing large scale structures such as 
voids or filament s. Note that there is evidence of a weak c lus- 
ter in the HDF (|Villumsen. Freudling fc da Costalll997l ). 

To model the source distribution from the HDF and 
HDF-S data, we first make two cuts to the catalogues. Ini- 
tially, we drop all objects with magnitudes greater than 
AB(8140) = 28.00. Next we select objects in the redshift 
interval 0.4 ^ z ^ 5.0. The lower bound is chosen because 
the lensing efficiency is negligible for low source redshifts. 
The upper redshift bound is chosen because the photomet- 
ric redshift determination is unreliable for z > 5.0 and most 
observed multiply imaged sources behind massive clusters 
have redshift z < 5.0. This yields 1688 objects when the 
catalogues are merged. The corresponding projected source 
number density is 215 arcmin -2 for the ensemble catalogue. 

We now bin the merged catalogue in magnitude from 



N(m,z) = / exp[— a (m — b) a ] exp[— c [z 
x exp[— g (m — b) u (z — d) v ], 



(9) 



where the exponents in the exponential terms, (s,t,u,v), 
can take an even integer value. The constants a, b, c, d, and 
/ are determined in the fitting procedure. We simply use the 
lowest \ 2 fit to construct a probability distribution function 
in the simulations. This is (s,t,u,v) = (4,2,8,2), with the 
other parameters given by: a = 1.607 x 10~ 4 , b = 30.70, 
c = 9.97x 10" 2 , d = 0.507, g = 2.197x 10" 7 and / = 18.462. 

Notably, the \ 2 °f the above fit is approximately one 
half of the \ 2 that is obtained when the binned magnitude 
and redshift data are decoupled and fitted on an individ- 
ual basis, i.e. N(z) and N(m). It is therefore clear that to 
avoid a systematic bias in the source generation, a joint fit 
in magnitude and redshift space is required. 



4-1.1 Comparison to Abell 1689 

Deep ACS imaging of Abell 1689 has allowed identifi- 
cation of_jnore__tlnin_30^familtes_of_jnu^^ imag es sys- 
tems (|Broadhurst et al.l liooBI ; iLimousin et all 120071 ). It is 
clear that Abell 1689 is the ideal candidate to compare 
with our source distri bution. We base our analysis on 
iBroadhurst et al. I (|2005l ). however our conclusions are un- 
changed if we use the analy sis of ILimousin et all ||2007| ). 
With the lens model as a guide IBroadhurst et al. ( 20051 ) have 
identified multiple images to / = 29.82. This means images 
with magnitudes comparable to and less than the magni- 
tude limit of our HDF catalogues, AB(8140) = 28.00, are 
detected in the ACS observations. This means we can com- 
pare our source distribution to the properties of Abell 1689 
directly. 

In order to compare the source characteristics of Abell 
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Figure 2. Comparison of image properties using the lensed HDF derived distribution and the Abell 1689 data of lBroadhurst ~et al.l 
(2005). In each bin, the simulated data has a mean given by the solid horizontal line and the gray rectangle represent ±lo\ Left panel: 
Redshift distributions. The Abell 1689 data are plotted as solid diamonds. Right panel: Magnitude distributions. The solid squares are 
the Abell 1689 data. 



1689 with our source distribution, we need to lens t h is dis - 
tribution through a mass model. Bro adhurst et all (|2005l ) 
have published an Navarro, Frenk, & White (NFW) fit (see 
Section r4.3.1|l to approximate the true mass profile. Using a 
radially averaged (and therefore circular) mass profile, they 
provide the tangential and radial critical curves, computed 
at redshift z = 3.0. They also determine the NFW scale ra- 
dius r B = 473j^gg kpc and we adopt r B = 473 kpc. We then 
adjust an NFW fit until the two critical curve locations are 
reproduced. 

Since we are using a smooth profile, the properties of 
Abell 1689 we compare our source distribution with cannot 
depend on the small scale details of the mass distribution. 
There are two quantities which satisfy this requirement: 1) 
the total number of lensed sources and 2) their redshift dis- 
tribution. A test of the HDF constructed source distribution 
would be the ability to reproduce these two quantities for 
simulated lenses similar to Abell 1689. Note that we do ex- 
pect some systematic effects due to the lack of assumed sub- 
structure in the modeling. In particular, the fact that Abell 
1689 is non-circular, and is known to have images which are 
influenced by substructure, means any comparison between 
the actual data and our model (which is circular, and has 
no structure) will exhibit residuals. 

To compare our simulations with the Abell 1689 
sources, we take Z beat as the source redsh ift for each im- 
age from table 2 of iBroadhurst et al.l |2005). Also, because 
the HDF catalogues used to derive the source distribution 
have been cut at z = 5.0, we drop all sources with redshift 
z > 5.0. The Abell 1689 source redshifts are then binned 
in intervals of Az = 1.0 from z — 0.0 to z = 5.0. For 
proper comparison with the simulations, the bin counts are 
normalized by the total number of families in Abell 1689. 

We simulate the smooth Abell 1689 NFW profile for 500 
source configurations generated with our source distribution. 
We then lens the sources independently for each configura- 
tion, using the procedure explained in Section [5T] Multiply 



imaged sources are identified and binned identically as the 
Abell 1689 sources. 

The average number of multiply imaged sources in the 
500 simulations is N s i m — 16 ± 4 per simulation. The min- 
imum and maximum number of families are 7 and 28 re- 
spectively and of the 500 configurations 102 exhibit 20 or 
more families of images . To z = 5.0 there are 27 families 
in the IBroadhurst et al.l (|2005l ) analysis. Although the mean 
number for our generated images is lower than that of the 
27 strongly lensed sources in Abell 1689 by ~ 3cr, we do 
not consider this a serious problem. The total number of 
strongly lensed sources is within factor of 2 in 69 per cent 
of the 500 realizations. 

We now examine the redshift distribution of the lensed 
sources in the simulations. The configurations are combined 
to produce a distribution of redshift probabilities in each 
bin. This allows the mean and standard deviation to be 
computed, see Fig. [5] Our simulated redshift distribution lies 
within la of the Abell 1689 sources, except for the ^ z < 1 
bin for which there are no observed multiply imaged sources 
in Abell 1689. 

Finally we discuss the image properties themselves. 
Since the mass model does not include substructure, we 
do not expect to reproduce the source multiplicity distri- 
bution or the image magnitude distribution of Abell 1689 
accurately. For a circular NFW mass profile, lensed sources 
with multiplicity of at most 3 can only be produced, whereas 
higher multiplicities are observed in Abell 1689. However, as 
shown in Fig. [2] the image magnitude distribution as gen- 
erated through our source distribution does match that of 
Abell 1689, and the overall shape is reproduced. 

Both 1) the total number, and the 2) redshift distribu- 
tion of the simulated sources are consistent with the source 
characteristics of Abell 1689. Given the results of this com- 
parison, we conclude that the HDF derived source popula- 
tion distribution is a faithful representation of the actual 
source distribution behind massive lensing clusters. 
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4.2 Cluster redshift distribution 

In order to generate a set of massive, strong lensing clus- 
ters in our simulations, we need to derive redshifts from an 
observed cluster redshift distribut ion. To do so, we u se the 
MAssive Cluster Survey (MACS) labeling et al.ll200lh sam- 
ple as our template. We choose this sample since it is X-ray 
luminosity selected, which is a robust proxy for mass, and 
the HST follow-up of these clusters has de monstrated th at 
a large majority of these are strong lenses (Ebcling 2005). 

We use the binned MACS cluster redshifts (see fig. 9 
m lEbeling et ail |200l])) to assign N\ ens (z) to each redshift 
bin. The z = 0.000 — 0.025 bin is omitted in the analysis 
since the cross-section for multiple imaging by clusters in 
this redshift range is negligible. An exponential fit is per- 
formed and the result is N\ ena (z) = 136.07 exp[— 5.845 z\. 
The fit is then normalized, integrated and inverted to form 
a probability distribution. The final probability distribution 
returns a cluster lens redshift for use in the simulations in 
the range z = 0.025 — 0.600 and can be expressed analyt- 
ically as P lcns (x) = -0.171 log 10 [-0.83406(- 1.0359 + as)], 
where a; is a random number between and 1. 



the range 10 14 — 10 15 Mq uniformly. We then assign pro- 
file parameters to characterize the mass distribution. A con- 
centration is drawn randomly from the range c = 4 — 5 
l|Bullock et al.ll200ll ) . Strong lensing models using NFW fit- 
ting profiles find that the typical scale radius lies in the 
range r s — 200 — 400 kpc. From the total mass, c and r s , we 
compute the characteristic velocity dispersion of the NFW 
halo. Note that we ignore any scaling of po with cosmology 
for simplicity. Finally, with the halo structure determined, 
the ellipticity is then drawn randomly from e ^ 0.15. 



5 ANALYSIS OF THE SIMULATIONS 

In this section we discuss the details of our simulations. 
Firstly, the procedure used to generate the lensing clusters 
and multiple image catalogues are discussed, and we explain 
how the source magnitude and redshift distribution and the 
cluster redshift distribution are used in the simulations. We 
then discuss our mass modeling assumptions. Finally, our 
marginalization method used to construct the confidence re- 
gions for cosmological parameter recovery are presented. 



4.3 Cluster mass distribution 

4.3.1 The NFW mass profile 

The NFW density profile is derived from fits to clus- 
ter scale dark matter haloes that form in A^-body 
simulations of a cold dark m atter dominated Universe 
(INavarro. Frenk fe WhitJl997f ). This density profile is given 
by, 



p(r,z) = 



Po 



((r 2 + z 2 ) 1 / 2 /r s )(l + (r 2 + z 2 ) 1 / 2 /^) 2 



(10) 



where z is orientated along the line of sight and the radius 
r is orthogonal to the line of sight. There are two defining 
profile parameters for the NFW model, the characteristic 
density po and the scale radius of the profile r s . The corre- 
sponding characteristic velocity dispersion v is defined as, 



2 4 2 
v = -Gr s p . 



(11) 



In our simulations we will calculate the total projected mass 
along the line of sight within the radius R = cr s , where c 
is the concentration parameter. For R > r s this is given by, 



M(R) 



3nv 2 r B 
2G 



In 



+ ■ 



(12) 



where R = cr s (see iLimousiru Kneib fc Nataraianl l|2005l ) 
for further details on the profile). 

4-3.2 Masses of the lensing clusters 

Strong and weak lensing analysis of massive cluster lenses 
have determined that the total mass of such cl usters within 
~ 1 Mpc, lies in the range, 10 14 - 10 15 M W , jSmith et all 
120051 ; lLombardi et al.1 120051 ; iNataraian et all \200H ). This 



mass range is the starting point for our construction of sim- 
ulated smooth cluster lenses. 

Consider a cluster scale halo with an NFW profile. To 
generate the halo we first draw a total cluster mass from 



5.1 Cluster and lensed image generation 

To perform our analysis, we generate clusters that will pro- 
duce the requisite number of families of multiple images. The 
methods used to generate the cluster lenses are described in 
detail in Section [4.3.21 We choose all clusters that produce 
multiple images which can be detected via relative astrom- 
etry in the High Resolution Channel (HRC) of ACS. This 
means the effective Einstein radius _Re, which we compute 
at z = 6.0, must be less than 50". In addition to this up- 
per bound we also apply a lower bound of _Re ^ 10". This 
bound restricts the simulations to massive clusters that plau- 
sibly exhibit image configurations that are measurable with 
the HRC. We note there is no effect on the simulations if 
the maximum Einstein radius condition is dropped, as long 
as the same pixel resolution is maintained. This may be the 
case if observations were taken with the Wide Field Cam- 
era on HST and the images were processed with sub pixel 
dithering. 

Once an appropriate cluster has been generated a source 
catalogue is generated. This catalogue is drawn from the 
HDF distribution derived in Section l4~Tl A redshift and mag- 
nitude is assigned for each source using the distribution de- 
rived (numerically) from equation ((9|. Each source is also 
assigned a random position in the field of view 47?e x 4_Re- 

The number of sources in the catalogue for which mag- 
nitudes, redshifts and positions are generated, is dictated by 
the projected source number density of the ensemble HDF 
catalogue which is 215 arcmin" 2 and the field of view. The 
source catalogue is then lense d with the code LENSTOOL 
|Kneiblll993l ; |Jullo et al.ll2007h . The lensed catalogue is then 
sorted into single images and multiple image families. 

At this point we have a complete simulated data set 
that is equivalent to a known observed sample of ideal lens- 
ing clusters. Observationally, this can be accomplished effi- 
ciently through high resolution space-based images of known 
massive clusters or deep ground based observations of clus- 
ters. We note that the MACS cluster survey provides a 



catalogue of clusters selected through their X-ray luminos- 
ity. This cluster catalogue provides a useful starting point 
for observing programs to identify strongly lensing massive 
clusters. Such a systematic survey is currently underway 
|Ebelindl2005h on EST ACS. 

5.2 Image catalogue construction 

The image catalogues we use in our recovery of the cosmo- 
logical parameters assume EST quality spatial resolution 
and spectroscopic redshifts. We first examine the ideal case, 
where the image positions and redshifts are known exactly. 
In Section 8 we expand our simulations to examine the ef- 
fects of errors in the image positions and source redshifts. 
However, before considering possible sources of error, a best 
case scenario must be established when cluster lenses are 
stacked. 

In the simulations, a magnitude limit is specified TOii m , 
when constructing the image catalogue from the lensed 
sources and is taken to be, u nless otherwise stated, mn m = 
24.5 in the AB(8140) b and l|Golse. Kneib. fc Soucaill 12002] ; 
iBroadhurst et al.l l2005f) . This ensures that all images will 
have measurable spectroscopic redshifts. Therefore, all fam- 
ilies in the lensed image catalogue which do not have two 
images with magnitudes m < mn m are removed at this 
point. 

Applying the magnitude cut prior to the final selection 
of families in the image catalogue is a stringent requirement, 
since we only use images with m < mn m in a given family. In 
the case of observational data, although it may not be possi- 
ble to identify all members of a given image family spectro- 
scopically, it is often possible to use properties like the image 
parity, morphology, photometric redshifts, or a combination 
of these to make the identification. Fina lly, since core im- 
ages have been detected in lensing studies (|Winn et al ]|2004 
IBroadhurst et al.ll2005l ). we do not remove such images un- 
less they do not pass the magnitude cut. 

5.3 Mass modelling 

In the spirit of our simulations, where we are investigating 
the ideal case of stacking strong lensing clusters, we use the 
same form of the mass profile that was used to generate the 
image catalogue. The centre of the lens is fixed at the origin 
and the ellipticity and lens redshift are taken as known from 
the original cluster lens. In Section [7] we describe how the 
profile characteristics of the mass modelling influence the 
cosmology recovery. It is shown in Section [7] that for the 
NFW profile, the velocity dispersion and scale radius should 
be considered unknown in the cosmology recovery. 

We note that the prior of a known mass profile is not the 
case when using observational data, where the underlying 
mass model is unknown. A variety of mass model priors will 
need to be considered in practise, since the mass distribu- 
ti on can be complicat ed, see for example the Bullet Cluster 
of lBradac et al.l (|200fJ ). When observational data is analysed 
the cosmological parameter bias from the modelling of real- 
istic clusters must be understood. These additional aspects 
present significant computational challenges when stacking 
a large number of clusters and we defer their study for the 
time being. 
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Table 1. Recovery Properties on the !!m-wx plane. The bounds 
are the limits of the lcr confidence contours, with the input cos- 
mology (!!m,wx) = (0.27, —1.00) subtracted. A\ is the area en- 
closed by the la contour (see equation 1141 ). 



Families 


Clusters 


Qm Bounds 


wx Bounds 




2 


5 


(-0.27,0.37) 


(-1.00,0.15) 


0.1825 




10 


(-0.20,0.26) 


(-0.56,0.11) 


0.0855 




10 


(-0.10,0.14) 


(-0.21,0.08) 


0.0227 




100 


(-0.07,0.08) 


(-0.10,0.05) 


0.0061 


3 


5 


(-0.26,0.21) 


(-0.60,0.20) 


0.1660 




10 


(-0.14,0.14) 


(-0.29,0.13) 


0.0447 




40 


(-0.09,0.08) 


(-0.12,0.07) 


0.0085 




100 


(-0.05,0.05) 


(-0.06,0.04) 


0.0025 


5 


5 


(-0.11,0.11) 


(-0.17,0.08) 


0.0149 




10 


(-0.08,0.08) 


(-0.12,0.08) 


0.0094 




10 


(-0.04,0.04) 


(-0.06,0.04) 


0.0023 




100 


(-0.02,0.02) 


(-0.02,0.02) 


0.0004 



5.4 Marginalization and cosmology recovery 

The % 2 minimization we use is computed in the source plane 
for computational speed, since the lens equation is many to 
1 when mapping from the image plane to the source plane. 
Our measure is given by, 

x =£g^sr- (13) 

where /3y is the source position corresponding to image j 
from the i th family as mapped back through the lens equa- 
tion (equation Q~|), and {5f3 2 )ij is the standard deviation 
squared in the source plane. The mean source position, (/3)j, 
is defined as (f3)i = Y17=i fiijl ni an d i s the average over all 
source positions of images belonging to the i th family. There- 
fore (f3 i; j — {/3)i) is the difference between the observed image 
positions as viewed from the source plane. 

In the calculation of equation (|13p . {5f3 2 )ij is computed 
using the assumed error in the image plane of 56 — 0.1". Us- 
ing the known image position, we sample at eight image error 
points: {6 x ,6 y ±56), {6 x ±56,6 y ), {6 X ± 56/ sfa, 6 y ± 50/y/2) 
and {6 X =F S0/V2, 6 y ± 56/y/2), and map each point back to 
the source plane using the lens equation. We then calculate 
the difference between /3 i • and each of the mapped error 
points. The eight 5/3 are first squared before averaging, i.e. 
{5[3 2 )ij, because this over estimates the standard deviation 
squared relative to the reverse order. Therefore using this 
averaging method leads to a conservative estimate of the 
final x 2 an d thus the confidence intervals (CIs). 

The minimization routine, which uses a parabolic op- 
timization in combination with parameter space bounds, is 
designed to ensure the x 2 surface of each cluster is accurate 
to at least 5x 2 < 0.002. This accuracy is sought because we 
want the final likelihood accurate to < 1 per cent. After this 
is achieved we form the la to 4a CIs. 



8 Gilmore & Natarajan 



1.0 

0.8 
0.6 
0.4 
0.2 
0.0 




1.0 
0.8 
0.6 
0.4 
0.2 
0.0 





0.0 



0.2 0.4 0.6 0.8 1.0 



0.0 0.2 0.4 0.6 0.8 1.0 



0.0 0.2 0.4 0.6 0.8 1.0 






0.2 0.4 0.6 0.8 1.0 



0.0 0.2 0.4 0.6 0.8 1.0 



0.0 0.2 0.4 0.6 0.8 1.0 

--M 




0.0 0.2 0.4 0.6 0.8 1.0 



(1.0 



-0.5 



-1.0 



-1.5 



-2.0 




0.0 0.2 0.4 0.6 0.8 1.0 



0.0 



-0.5 



-1.0 



-1.5 



-2.0 



5 families 
5 clusters 




0.0 0.2 0.4 0.6 0.8 1.0 



-0.5 




-2.0 



0.(1 



-0.5 



-1.0 



-2.0 





3 families 


0.0 


i i i i i i 


5 families 




10 clusters 


-0.5 




10 clusters 






£ -1.0 










-1.5 










-2.0 















0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 

I\i 



Figure 3. Confidence contours on the Om-^A and !1m-wx planes for 2, 3 and 5 multiple image families (columns 1, 2 and 3 respectively). 
For each family number specification, 5 or 10 clusters are stacked. The input cosmology marked by the cross in all panels, is a flat ACDM 
universe with (Omi^X; w x) = (0.27,0.73,-1.00). In the top two rows wx = —1.00 and a flat universe is not assumed during the 
cosmology recovery. In the bottom two rows a flat cosmology is assumed with wx and S!m used in the cosmology recovery. In all panels, 
the recovery is undertaken with both the velocity dispersion and the scale radius allowed to vary. The ltr confidence region is the internal 
white region centred on the input cosmology. The 2a, 3cr and 4cr confidence regions are progressively darker. 
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6 RESULTS: RECOVERY OF INPUT 
COSMOLOGY 

We investigate the cosmology recovery properties by con- 
sidering cluster samples with different family numbers. To 
obtain the recovery properties, we proceed as outlined in 
Sections [4] and [5] At this stage the input parameters to the 
simulations are the number of clusters to stack, the total 
number of families in each cluster and the magnitude limits 
for the multiple images. The recovery properties are com- 
puted for each cluster, in the specified cosmological param- 
eter space. The constraints from each cluster in the sample 
are then stacked to obtain the final recovery properties. 

We consider observations of 5 and 10 clusters, with 2, 3 
and 5 multiple image families. The magnitude limit, mn m , 
is AB(8140) = 24.5, 24.5 and 26.0 for the 2, 3 and 5 multi- 
ple image families respectively and the recovery properties 
of the input cosmology are computed in the Om-^a and 
!!m-wx planes. The results of stacking the first 5 generated 
clusters and then all 10 are shown in Fig.[3]and characterized 
in Table Q] for the !Jm-wx plane. 

In the f^M-fiA plane we observe that the parameter con- 
straints obtained by stacking lensing clusters is insensitive 
to the value of the cosmological constant when marginaliz- 
ing over Qm- The Qm constraints are also not stringent, with 
Qm constrained to lie between approximately 0.1 and 0.50 
for the 5 family 10 cluster case. 

These properties of the cosmology recovery in the Qm- 
Qa are expected. By examining Fig. [T] one can see that the 
contours passing through the input cosmology are degener- 
ate for different source redshift configurations. This means 



that the constraints derived from the family ratio will be 
tightest around the contour defined by the input cosmology. 

However, in the Qm-wx plane (lower two rows of Fig. [3j 
the family ratio is sensitive to the dark energy equation 
of state. In Fig. [T] we see that the S contours passing 
through the input cosmology are quasi-orthogonal for dif- 
ferent source redshift configurations. Therefore clusters with 
distinct source redshift planes will have different constraint 
directions in the S7m-wx plane. This is also the case for dif- 
ferent lens redshifts given the same two source redshifts how- 
ever the effect less pronounced. When stacking the cluster 
sample the orthogonality of different family ratios leads to 
closed contours in the S)m-wx plane, for the 5 family 10 
cluster simulation. 

We then examine the constraints obtained when 40 and 
100 clusters are stacked for 2, 3 and 5 multiple image fami- 
lies. The magnitude limits are the same as for the 5 and 10 
cluster runs but only the SIm-wx plane is considered. The 
results are shown in Fig. [4] and Table [1] Note that 40 clus- 
ters with 5 families and 100 clusters with only 3 families, 
produce almost identical CIs. 



7 DEPENDENCE ON LENS PARAMETERS 

The role of the lens parameters in the optimal retrieval of 
the input dark energy equation of state is examined in this 
section. We address the issue of the \ 2 dependence on the 
marginalization over combinations of the dominant lens pa- 
rameters. The marginalization for each parameter set is com- 
puted on the f^M-wx plane in a flat Universe, for an input 
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cosmology (Om,Wx) = (0.27,-1.00). We also consider the 
constraints derived from a full multi-dimensional marginal- 
ization and compare to the corresponding simple \ 2 mini- 
mization. 

7.1 Marginalization over lens parameters 

A massive NFW halo is considered in this section and sub- 
structure is omitted from the computations for the purpose 
of simplicity. The NFW halo is fixed at the origin and the 
parametrization we us e is typical of massive cluster lenses 
l|Nataraian et al1l2007t ) and the lens parameters we use are: 
v = 1600 km s _1 , r s = 200 kpc and e = 0.10, with the lens 
redshift z\ = 0.3. Two sources are lensed and they are lo- 
cated at (x i, yi) = (-1.00, 0.10) and (3:2,2/2) = (1.20,4.80) 
respectively, with redshifts z s \ = 1.0 and z S 2 = 3.0. Sources 
1 and 2 have multiplicities of 5 and 3 respectively. 

7.1.1 Case with known lens parameters 

First we consider the case where one has complete knowledge 
of the gravitational lens, i.e. all lens parameters are known 
exactly. The result of the marginalization with all lens pa- 
rameters fixed is shown in the upper left panel of Fig. [5] for 
the NFW lensing cluster. Note that the confidence region 
for the NFW profile is non-zero. In fact, the existence of a 
non-zero confidence region is a general feature of the lensing 
systems in our simulations. 

Based on Fig. [S] it is clear that the intrinsic confidence 
region of a lensing system represents a fundamental limit on 
the determination of the cosmological parameters from a sin- 
gle configuration. Even in the optimal scenario, i.e. complete 
knowledge of the lens system, one must consider the confi- 
dence region and not the best fitting cosmology in a rigorous 
statistical treatment. The best-fitt ing cosmological param- 
eters, as for example employed in iDalal. Hennawi fc Bodei 
(2005), cannot be used to characterize the cosmological pa- 
rameter recovery, since they do not take into account the 
limitation of the intrinsic confidence region. Noting this 
property, we have used CIs exclusively in this work. 

7.1.2 Mass profile parameters 

In practice, one does not have prior information on all the 
parameters defining the mass profile of an observed lensing 
cluster. For the NFW lensing profile, there are three com- 
binations of the lens parameters to consider marginalizing 
over when recovering the cosmology. These are: the veloc- 
ity dispersion v, the scale radius r a and both v and r s . The 
results of each marginalization are given in Fig. 

For the NFW profile, the shape of the confidence con- 
tours changes for the different marginalization parameters. 
For the case of marginalization over r s alone, in which v is 
held fixed, closed confidence contours are obtained in the 
S7m-wx plane. This is distinct from the contours derived for 
the marginalization over v alone, which yields a confidence 
region which follows the contours of constant cosmological 
family ratio 5. In Fig. [5] the simultaneous marginalization of 
r s and v expands the 4a confidence region asymmetrically. 
In addition, note that all other confidence contours qualita- 
tively follow the contours for the family ratio. Therefore it is 
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Figure 5. Marginalization contours over different NFW param- 
eters. The la- to 4<r contours in the Qm~ w x plane are shown, and 
the cross marks the input cosmology (S!m,wx) = (0.27,-1.00). 
The parameter(s) in the upper right corner of each panel specify 
the marginalization variable(s). 

obvious that optimization over both the characteristic den- 
sity and scale radius is required when recovering cosmology 
constraints from the NFW profile. 

An additional comment on Fig. [5] is warranted with re- 
gard to the cosmology recovery. It is obvious from the panels 
in Fig.[5]that when only r s or only v are allowed to vary, that 
tighter cosmological constraints are obtained, compared to 
the case where both r s and v are marginalized over simulta- 
neously (bottom right panel of Fig. [5}. It is therefore natural 
to ask, what are the cosmology recovery properties in the 
cases where only r s or v are marginalized for stacked clus- 
ters? However, when considering actual cluster lenses, if we 
were to marginalize over say only r s , we would have to fix the 
velocity dispersion. This would require an additional prior, 
and in actual cluster lenses this prior is not typically avail- 
able. It may be possible to constrain the total cluster mass 
through X-ray observations, and in this case, one would use 
equation (|12[l to apply the priors to both the velocity disper- 
sion and the scale radius. However, in general, when using an 
NFW profile to model an actual strong lensing cluster , both 
r s and v are not known a priori dComerford et al.ll2006h . The 
reader should therefore note that although stronger cosmo- 
logical constraints could be obtained if r s or v were known a 
priori, this is not the case in actual cluster models because 
separate priors on r s or v cannot be found. To this end, we 
have only presented the marginalization over both r s and v 
in Section [6] 

7.2 Full marginalization versus x 2 minimization 

Up to this point, we have considered a simple minimization 
routine of the lens parameters with respect to the \ 2 estima- 
tor. We have done this over the cosmological plane and for 
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different lens parameter combinations for the NFW profile. 
However, we now examine the validity of this assumption 
for the NFW lensing configuration in Section \7. II 

We now perform a full multi-dimensional marginaliza- 
tion over the lens parameters, r s and v, for the NFW lens- 
ing configuration. Our adopted method is as follows: first, we 
find the optimized lens parameters which gives the minimum 
Xmin a t each cosmology point. From the optimized lens pa- 
rameters, we then compute a x 2 surface extending to \ 2 = 
(Xmin + 50). The parameter spacing on the \ 2 surface is 
0.2 arcsec for r B and 0.2 km s _1 for v. The \ 2 surface is then 
converted into a likelihood surface using, C = exp(— X 2 /2), 
which is integrated to obtain the marginalized likelihood at 
a given cosmology input. This likelihood is then converted 
to a x 2 surface, using x 2 = — 2 ln(£/£ max ). 

To quantify the difference between the multi- 
dimensional marginalization and the simple optimization 
method, we consider an integral which characterizes the en- 
closed area, 



An — 



/ ■ 

JR{n) 



(14) 



where n defines the confidence contour, and R(n) the region 
enclosed by the given contour. 

For the multi-dimensional marginalization over r s and 
v, we have approximately, A\ = 0.267, A2 = 0.462, ^3 = 
0.680 and A4 — 0.889. Whereas for the simple x 2 minimiza- 
tion Ai = 0.273, A 2 = 0.464, A 3 = 0.683 and A A = 0.890. 
The difference between the enclosed area is largest for A\, 
which is ~ 2 per cent. Importantly, the la contour from the 
simple optimization method encloses an area larger than 
the full multi-dimensional marginalization. Both methods 
recover a x 2 — at the input cosmology. 

Given the maximum ~ 2 per cent difference between 
the enclosed areas from the multi-dimensional marginaliza- 
tion and the simple x 2 -°ptimization, we are justified in us- 
ing the simple optimization method. With sufficient multiple 
images to constrain the cluster potential, one expects only 
one global minimum in the x 2 optimization surface. This is 
the case for the NFW lensing configuration considered here, 
and for similar lensing configurations one expects the op- 
timization method to accurately reproduce the confidence 
contours derived from a full multi-dimensional marginaliza- 
tion. 



8 STUDY OF SURVEY CHARACTERISTICS 

We have also simulated smooth cluster lenses that contain 
observational noise. The noise sources we have considered 
are positional errors and redshift error. These two sources 
are the two major observational parameters which will deter- 
mine the type of survey required for the application of CSL. 
Our simulations are restricted to the 5 families 10 cluster 
case (see bottom right panel in Fig. [3} in the J7m-wx plane. 

Our objective is to examine the levels of noise which 
do not lead to relaxation of the CIs, when compared to our 
simulations of the ideal case. We use this approach since 
this will represent a tolerance threshold where parameter 
constraints will weaken. 

First we consider redshift accuracy. We have re-run 
the 5 families 10 cluster simulation with gaussian red- 



shift noise with the following standard deviations Az = 
0.003,0.005,0.01, and 0.1. When adding the noise to the 
source redshift, we draw from the appropriate gaussian and 
add the resultant value to the known redshift. Note that 
all multiple images of a given family receive the same noise 
contribution to their redshift. Our results are as follows. We 
find that photometric redshift accuracy, i.e. Az — 0.1, leads 
to significant relaxing of the CIs and in some cases does 
not recover the input cosmology correctly. Of the remaining 
standard deviation simulations, all recover the input cosmol- 
ogy with increasingly tighter constraints. For Az — 0.003, 
there is negligible confidence contour relaxation when com- 
pared to the ideal simulation. 

Next we move onto image positions. Following a similar 
procedure to the redshifts, we find that a positional accu- 
racy of ±2 pixels of the HST ACS is required to recover 
the confidence region in the SIm-wx plane. This means that 
space based imaging of the strong lensing clusters will be 
required. We have also considered image position error in 
combination with redshift error, and find that the combi- 
nation Az — 0.003 and ±2 pixels in ACS does not lead to 
relaxation of the confidence contours compared to the ideal 
case. This means that space based imaging of strong cluster 
lenses, followed up by spectroscopic redshift determination 
of the multiple images, is the appropriate survey method to 
pursue for the application of CSL to cosmography. 



9 LIMITATIONS OF OUR CURRENT 
ANALYSIS 

As with any other method used to constrain cosmological 
parameters, CSL has a number of sources of error one must 
consider. Our goal has been to examine observational strate- 
gies under the most optimal conditions in this first detailed 
study. Our results demonstrate that by stacking lensing clus- 
ters with several multiple image families, confidence regions 
for the input cosmological parameters can be recovered, as 
shown in Figs [3] and [J] 



In earlier work, iDalal. Hennawi fc Bodd l|2005h raised 
the issue of parameter bias introduced by the unknown clus- 
ter mass profile, and the presence of substructure. They con- 
cluded from extracting the best-fitting cosmological param- 
eters (rather than confidence regions) that this technique 
produces > 100 per cent errors for a single cluster. We 
are not disputing the importance of correctly modelling the 
mass profile and substructure. However, in our simulations 
we have focused on the recovery properties of the ideal sce- 
nario, and the errors cosmography CSL can tolerate in an 

obs ervational data set. 

IDalal. Hennawi fc Bodd l|2005l ) has also investigated the 
recovery bias including line of sight density fluctuations us- 
ing a single numerical cluster. It should also be pointed out 
that any parameter recovery bias resulting from substruc- 
ture a long the line of sight (fig. 4 of lDalal. Hennawi fc Bode 
will be minimized if uncorrelated lines of sight are 
used when constructing the cosmological constraints. Also, 
bias due to line of sight substructure could possibly be ac- 
counted for, by performing redshift surveys behind cluster 
lenses and modeling any detected structure along the line 
of sight in multiple lens planes. Removing the systematic 
parameter bias due to substructure along the line of sight is 
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an important topic that requires further investigation and 
is left for future work. 

Since we have only used positional accuracy in our cos- 
mological recovery, we should ask: What additional \ 2 mea- 
sures can we add to our current recovery measure of equation 
(|13p that will further constrain the input cosmology? The 
first obvious candidate is a higher order shape estimators for 
multiple images, since we have ignored all shapes and simply 
used a point approximation in our simulations. In real clus- 
ters multiply imaged sources will have their shapes distorted 
as they are strongly lensed and this will give additional infor- 
mation about t he len sing mass distribution, see for example 
iLeonard et al.l (|2007T ) . Another improvement may come from 
the possibility of using other observational maps of clusters 
to infer additional properties of the cluster. For example, 
a total mass constraint can be obtained from X-ray cluster 
observations. 



10 CONCLUSIONS 

In this paper, we have examined the recovery of input cos- 
mological parameters using the cosmological family ratio for 
several sets of multiple images in stacked simulated smooth 
lensing clusters. For a given lensing configuration, the con- 
straints on the JIm-wx plane are tighter than those on the 
Om-^a plane. We find that the recovery of input parameters 
depends both on the number of multiple image families with 
known redshifts and the number of clusters stacked. 

As shown in Figs [3] and [3] the confidence contours 
tighten as we increase the number of families from 2-5, and 
as the total number of stacked clusters is increased from 1 
- 100. For example, the recovery is comparable in terms of 
area enclosed within the la - 4a contours for the case with 
5 families using 40 clusters and 3 families using 100 clusters. 
As we demonstrate, both strategies - stacking fewer clusters 
with a larger number of image families or a larger number of 
clusters with fewer families lead to comparable constraints. 

Currently there are about 40 clusters with high resolu- 
tion HST ACS images, of which a handful (about 7 - 10) 
have more than 3 families of multiple images with measured 
redshifts. Proposed future space missions, the James Webb 
Space Telescope and JDEM will enlarge the existing lensing 
cluster samples by a factor of 20 or more, while simultane- 
ously detecting a larger number of families on average for 
each cluster through very deep imaging. Determining spec- 
troscopic redshifts for multiple image families is likely the 
limiting step in terms of the observational effort for the CSL 
method. 

If 10 clusters with 3 multiply imaged sources are 
stacked, the dark energy equation of state in a flat universe 
can be constrained to —1.6 < wx < —0.8 at 2a where Qm 
has been marginalized over, see Fig. [3] This is comparable 
with current methods (|Komatsu et all 2009). Therefore to 
be competitive with current methods, given that actual mas- 
sive cluster lenses will present sub-optimal recovery condi- 
tions relative to those considered here, an initial application 
of CSL to actual lensing clusters will require at least 10 clus- 
ters with 3 or more multiply imaged sources with Abell 1689 
quality HST ACS images. Such constraints may be obtained 
only if the source redshifts are determined spectroscopically. 
These factors indicates that initial observations should focus 



on massive clusters that exhibit large lensing cross sections. 
Once these strong lensing clusters are imagined to Abell 
1689 quality with the HST, those with 3 or more multiply 
imaged sources should be given priority when determining 
which lensing clusters to follow up spectroscopically. 
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